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ABSTRACT Some features of hydro- and thermodynamics, as apphed to atmospheres 
and to steUar structures, are puzzhng: 1. The suggestion, first made by Laplace, that our 
atmosphere has an adiabatic temperature distribution, is confirmed for the lower layers, 
but the reason why it should be so is difficult to understand. 2. The standard treatment 
of relativistic thermodynamics does not allow for a systematic treatment of mixtures, such 
as the mixture of a perfect gas with radiation. 3. The concept of mass in applications of 
general relativity to stellar structures is less than completely satisfactory. 4. Arguments 
in which a concept of energy plays a role, in the context of hydro-thermo dynamical sys- 
tems and gravitation, are not always convincing. It was proposed that a formulation of 
thermodynamics as an action principle may be a suitable approach to adopt for a new 
investigation of these matters. 

In this second article of a series we propose to base criteria of stability on the hamil- 
tonian functional that is provided by the variational principle, to replace the reliance that 
has often been placed on ad hoc definitions of the "energy". We introduce a new virial 
principle that is formulated entirely within the Eulerian description of hydrodynamics, 
which allows a simpler derivation of a well known stability criterion for polytropic stellar 
configurations. Boundary conditions are based entirely on mass conservation. 

The new approach is tested on isothermal and polytropic atmospheres and then used 
to initiate a new study of stars. Traditional results for polytropic, spherical configurations 
are confirmed, but our study gives new insight and new results for the case that radiation 
is taken into account. 

PACS Keywords: Atmosphere, astrophysics, action principle, virial theorem. 
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I. Introduction 



Stability of a dynamical system that consists of a gas held together by mutual gravita- 
tional forces can be studied by comparing the energy of a presumptive bound configuration 
with that of an infinitely diffuse state. But a necessary condition, for such comparison to 
carry any weight, is that a framework exists within which a well defined hamiltonian func- 
tional of the dynamical variables can be provided. 

The literature contains stability studies where this condition, if not clearly stated, 
was nevertheless satisfied. But there are numerous instances, even in textbooks, where an 
expression for the energy is written down without any supportive dynamical framework. 

In these papers we insist on basing all investigations on a variational principle. The 
dynamical equations are in most cases the same as those of previous studies, but in the 
case that the effect of radiation is included we are breaking new ground. 

The paper begins with a very brief summary of the first paper of this series (Fronsdal 
08). A well known variational principle for irrotational hydrodynamics is expanded by 
elevating the temperature to the role of an independent dynamical variable. After a brief 
summary of standard atmospheric models we introduce a new "virial" theorem. It is an 
analogue of the classical virial theorem, in which the particle coordinate operator r ■ p is 
replaced by p$, where p is a mass density and $ is the velocity potential (the canonical 
conjugate of p). Using this new virial theorem in a manner that recalls Jacobi's's original 
use of the classical analogue we are led to the well known stability criterion (n < 3) for 
polytropic atmospheres. No appeal is made to kinetic theory and no approximations are 
made. 

Section 2 investigates the stability of several atmospheric models. The boundary 
conditions, that are central to the question of stability, are deduced in each case from the 
conservation of mass. This is a reminder of the fact that the usual approach to relativistic 
thermodynamics abandons the continuity equation of classical hydrodynamics, a defect 
that can easily be repaired (Fronsdal 2007) . 

We study the stability of several atmospheric models by perturbation theory. The 
important role of mass conservation in the choice of initial conditions is emphasized. 

This paper is limited to the non relativistic approximation. The relativistic case will 
be examined in another paper of this series. 



1.1. Hydrodynamics 

Basic hydrodynamics deals with a density field p and a velocity field v over H^, subject 
to two fundamental equations, the equation of continuity. 



p + div(p^*) =0, p: 



dp 



(1.1) 



and the hydrodynamical equation (Bernoulli 1738) 




D 



(1.2) 



This involves another field, the scalar field p, interpreted as the local pressure. 
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We assume that the velocity field can be represented as the gradient of a scalar field, 

V — — grad $. (1.3) 

In this case the hydrodynamical condition is reduced to 

grad p — p grad ($ — jl)- (1.4) 

To complete this system one needs a relation between the fields p and p. 
A local functional V\p\ is defined up to an additive linear term by 

p = pV' - V, V' := dVjdp. (1.5) 

Then dp = p dV and the equation (1.4) becomes, if p 0, 

V' = + A, A constant. (1.6) 

The two fundamental equations that characterize irrotational hydrodynamics are 

p + div(pv) = 0, p:= dp/dt, 
^-v'^/2 + X = dV/dp, 

together with the defining equations v = —grad $, p :— pV — V. These equations are the 
Euler-Lagrange equations associated with the action (Fetter and Walecka 1980) 

A[p,^ = J dtd^x C, C = p{^-if/2 + X)-V[p]. (1.8) 

The value of this last circumstance in the present context lies in the fact that it gives us a 
valid concept of a total energy functional. 



1.2. The mass 

The conserved density p will be taken to have the interpretation of mass density, and 
the total mass is the constant of the motion 




Since the total mass is a constant of the motion it is natural to fix it in advance and 
to vary the action subject to the constraint J d^x p{x) = M. The parameter A takes on 
the role of a Lagrange multiplier and the action takes the form 

A = J d^x(^p{^-v'^/2-(l))-V^ +X(^J d^xp-My (1.9) 

The gravitational potential (j) is included. 

The conservation of mass has important implications for boundary conditions. 
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1.3. Isothermal and polytropic atmospheres 

An ideal gas at equilibrium, with constant temperature, obeys the gas law 

p/p = nr. (1.10) 

Pressure and density are in cgs units and /iTZ = .8314 x 10® erg/K, where ji is the atomic 
weight. 

In hydrodynamics, the isothermal atmosphere can be given a lagrangian treatment by 
taking 

V = nTp\ogp. (1.11) 

We suppose that the gas is confined to the section zq < z < zq -\- h oi vertical cylinder 
with base area A and expect the density to fall off at higher altitudes. The action density, 
for a gas at constant temperature T in a constant gravitational field = gz, g constant, is 

C[^,p\=p{(^-^^'^/2-gz + \)-nTp\ogp. (1.12) 

We may consider this an isolated system with fixed mass. 

At equilibrium $ = 0,i; = 0, p = and the equation of motion is V' — 7^T(1 + log p) = 
A — gz, hence 

p{x, y, z) = e-l+V7^Tg-gz/7^T^ ^ ^ ^^e-i+^/^'^(l - e-^VRT) ^-^^o/tzt 

g 

and after elimination of A 

gM e-g(^-^o)/''^'^ gM e-e(^-^o)/^'^ 

^ " ATZT 1 - e-gh/RT ' ^ = ^ i _ g-gVRT • (l-^^) 

(If the atmosphere extends to infinity the pressure at the bottom is gM/A, as it shoul be.) 
A difficulty with this model is that the specific internal energy is TZT log p, which is very 
different from that of an ideal gas. 

The isothermal atmosphere is usually abandoned in favor of the polytropic atmo- 
sphere. A polytropic gas can be described by the lagrangian (1.9), with 

V = ap'^ , a, 7 constant. 

Variation with respect to p gives 

a 1 
P^-P^: -=7-1. 
n n 

The temperature does not appear explicitly but is taken to be determined by the gas law, 
Eq.(l.lO). At mechanical equilibrium v = 0,p = and X — gz = a7p^/"', hence 

a7 
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Since the density must be positive one does not fix the volume but assumes that the 
atmosphere ends at the point zi = X/g. Then 



This fixes h and thus zi and A. If the atmosphere is an ideal gas then the temperature 
varies with altitude according to 

7^T = p/p= V^ = ^^ (1.14) 
n n + 1 

For air, with atomic weight 29, TZ = 2.87 x lO^erg/gK and n = 2.5. At sea level, 
g = 980cm/sec^, the density is p = 1.2 x 10~^g/cm^, the pressure p = 1.013 x lO^dyn/cm?. 
Thus p/p = .844 X lO^cm^/sec^, T = Tq = 294K, = 3.014 x lO^cm 30A;m. and the 
dry lapse rate at low altitudes is —T' = 294/^1 = 9.75K/km. 

The specific internal energy of this model is np/p — nTZT, as it should be for an 
ideal gas. However, the temperature is not treated as an independent dynamical variable, 
instead it is fixed by the constraint (1.10). The theory developed so far is partially "on 
sheU". 

A variational , principle is incomplete without the specification of boundary conditions. 
The only boundary condition that we shall apply systematically is the conservation of mass, 
which controls the material flow at the boundary. 



II. Thermodynamics 

2.1. The adiabatic lagrangian 

The action principle can be extended so as to treat the temperature as an independent 
dynamical variable. 

Two kinds of additions can be made to the lagrangian (1.9) without spoiling the 
equations of motion that are essential to hydrodynamics. To preserve the equations of 
motion as well as the correct expression for the internal energy, we have been led to the 
lagrangian for an ideal gas, 

p, T] = p($ - v'/2 -cp + X)-nTp log f + (2.1) 

with k := p/T"^- The last term is the contribution of radiation, assumed as is usual to 
be that of a black body; the constant a is the Stefan-Boltzmann constant, a = 7.56 x 
10-^^ ergs /cm^K^. 

Variation of the velocity potential $ gives the continuity equation as before. Variation 
of the density p leads to 

^ -v^/2- + X = nT(l+\og^). (2.2) 

ko 
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Taking the gradient one gets 



- pgrad (j) = -grad p, (2.3) 
where D/Dt — d/dt + \/ ■ v is the total derivative time and 

p = nTp+^T\ (2.4) 

Variation with respect to T gives 

n(n-log-)p+-T' = 0, (2.5) 

With the neglect of the radiation term, or if n = 3, this reduces to the polytropic equation 
of change, k = constant. 

For conditions as in the atmosphere of the Earth, log(/c//co) — n ~ 10"^*^. When the 
emphasis is not on applications we shall use units such that ko = 1. 

The on shell hamiltonian density, obtained from (2.1) with the help of (2.5), takes the 

form 

h = pv'^/2 + p(l) + u, u = nnTp + aT'^. (2.6) 

We suggest that, when n 7^ 3, using the lagrangian (2.1) is preferable to the usual 
assumption /? := /Jgas/ptot = constant, which is true only when n — 3. 

We have an action principle, with dynamical variables p and T , that reproduces all of 
the equations that characterize the equilibrium configurations, including the contribution 
of the photon gas to the energy and the pressure, as well as the standard, hydrodynamical 
relations of an ideal gas. 



2.2. Virial theorem 

Both (2.4) and (2.6) are usually derived from considerations outside the proper do- 
main of thermodynamics. We prefer an axiomatic foundation of thermodynamics that is 
complete in the sense that it does not need other input. As an example let us discuss the 
use of the virial theorem to make certain predictions concerning stability. 

The virial theorem was introduced into the present context by Kelvin. It is based on 
the scaling properties of the hamiltonian of a system of particles. li H = K + V, kinetic 
energy plus potential energy, then the lagrangian is K — V and the equations of motion 
imply that, up to a time derivative, 

^if-p) ^2K -f- vV, 

where p = mv, K is the kinetic energy ^/2 and V is the potential. If the potential is 
homogeneous of degree n then 

d 

— {r-p) = 2K -nV, 
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If the particle goes through a cycle then the average of this quantity over the cycle is 
zero. In the case examined by Kelvin the potential is homogeneous of degree -1, so that, 
in the case of periodic motion , when average is taken over a period, V = —2K. This is 
used to prove stability, since the total energy is E — K + V — —K is negative. Actually, 
what is proved is that, if the motion is periodic then it is bounded. 

The application to continuous systems, developed by Clausius (1851), is much more 
difficul. According to Chandrasekhar (1938) (pp. 49-51), who also quotes Poincare, the in- 
ternal energy is the kinetic energy associated with the microscopic motion of the molecules. 
It is assumed, usually without discussion, that the presence of gravitational forces do not 
affect the internal energy, and that the total energy is obtained by simply adding the grav- 
itational potential energy to it. In the present approach we are led to the same conclusion, 
the energy is identified with the hamiltonian. But the theory does not authorize drawing 
on kinetic theory to formulate a virial theorem. In fact, some attempts to do so have been 
successful only after making various assumption and/or approximations. See for example 
CoUins (2007). 

There is; however, a virial theorem associated with a lagrangian of the type (2.1), that 
we abbreviate as 

£ = p($ - ^-2/2) - V. 

(The potential V includes the gravitational field and the lagrange multiplier.) Variation 
of $ and of p gives the equations of motion 

^ = v'^/2 + {dV/dp), p = -div(p^T), 

In accord with the boundary condition discussed at the end of Section 1.3 (the conservation 
of mass), we drop the boundary term — J ^pv ■ da to get 

/.3.|(^) = /d».(pf-p.V2). 
If the system goes through a cycle then the average of this quantity over the cycle is zero, 

J dx pv^/2 = J dx p^. (2.7) 

In the case of the lagrangian (2.1), but neglecting the radiation term, we obtain, 

J dx (pv'^/2 + p(0 - A) + 7^Tp(l + log k)^ = 0. (2.8) 

To draw conclusions from this result we adopt the strategy of Jacobi (1889) , noting 
that the total energy is 



E = J d^x(^p{v^/2 + (j)/2) + nTplogky 
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By replacing the gravitational potential by </)/2 we are specializing to the case that the 
interaction is between the particles, not with a fixed source. In this case the interaction 
energy is bilinear in the density and variation of p introduces a factor of 2. The equations 
of motion are always the same but the factor 1/2 is needed in the action. 
Combining the last two equations we obtain 



7^= [ d^x(-pv^/2 + {- -l)p(f)/2 + Xp). 
J \n n ) 



The last two terms can be evaluated. The calculation was done first by Ritter (1870); it 
is in many textbooks, and it is repeated later in this paper. The first result. 



y d^xXp = XM = X(f){R) = 



R 



is obtained with the help of the boundary condition at the surface and elementary prop- 
erties of the gravitational interaction. The derivation of the following formula, 

i'.,^l2 = E, = ——. 

makes use the known degree of homogeneity (7) of the heat term with repect to the density. 
It is important that the equation of motion was not used. The sum is 



j d^x({^-l)p(t>/2 + Xp^ 



(n-3)(n + 1)^7^2 
(5 — n)n R 



This shows that, when n > 3, the potential energy near the point of equilibrium is positive, 
which allows the gas to diffuse; the configuration is unstable. If n = 3 there is a stable 
solution with zero energy, and another solution for which the energy is not defined. It just 
walks away, like a free particle or a critically damped harmonic oscillator. 



III. Stability of some atmospheres 

In all the examples of fluctuations around a static solution that follow we consider 
only adiabatic fluctuations; in the isothermal case the temperature is flxed and in the 
polytropic case the fluctuations of T are governed by 5p/ p = n5T/T, by fixing k = p/T'^. 
One has to ask if this limitation to adiabatic fluctuations is fully justifled. To study this 
question it is necessary to introduce entropy and to consider the heat equation, as we have 
attempted to do in the flrst paper (Fronsdal 2009). As pointed out by Emden (1907), it is 
necessary to assume that the heat flow has zero divergence, and this is supported by the 
fact that the polytropic atmosphere has constant lapse rate. 
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3.1. The isothermal column 

The equations of motion (1.7), with V = TZTplogp (see Eq.(1.12)) are 

p + (pv)' = 0, v+ + 7^T(1 + logp))' = 0. 

The prime indicates derivation with respect to z. We consider the space that is tangent 
to a static solution with density po- Setting p = pQ + Sp we have the foUowing equations 
for the perturbation Sp, 

6p+{povy, v + TZTidp/poY ^0, 

Thus 

Sp = 7iT{poa'y, a:=Sp/po. 
For a harmonic mode with frequency u>, 

-io^Sp = RT{poa'y, (3.1) 

and ^ 

"■^T J '^^ P^^^ ~ J '^{PoCi'y dz — pQaa' —J poO.''^dz. 

The configuration is stable if this implies that o;^ > 0, which will be the case if the 
boundary term vanishes. To justify any choice of boundary conditions we have only the 
conservation of mass, J Spdz — 0. This ensures that Sp fall off at infinity and we are left 
with -5p(0)a'(0). 

We shaU show that a'(0) = 0. Eq.(3.1) tells us that 

= (Po/Po)"' + oi" = + 
This is a linear differential equation with constant coefficients, with general solution 



a(z) = Ae^+^ + Be^-^, k± = ± \l (^—) 
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27^T V '27^T' 7^T 
Since, up to an irrelevant constant factor, po = 6xp(— ^2;/7?.T), 



Sp = poa = Ae^+^ + Be^-^, a± = — ^ ± \ {^-)^ 



^2 



27^T V 27^T' 7^T 

In order that SM = Aj Spdz vanish we need for the coefficients A, B to be both non zero, 
and convergence of the integral then requires that a;^ > 0. In this case 

r,. /"r , A B A B 

SM = / Spdz = = — + — , 

J a+ a_ k- k+ 
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the vanishing of which requires that a'(Q) = 0. Therefore, not only is the condition 
> verified; it is also confirmed that the boundary condition a'{0) = is the only one 
possible. We have seen that this choice of boundary conditions is the one that ensures the 
conservation of mass. 

The fluctuation 5p does not vanish at the lower end, but the velocity does, as is 
natural. In fact, we can reach the same conclusion more easily as follows. Conservation of 
the current implies that 

CXD 

pv = 0. 



Since, as we showed, the upper limit cannot give a contribution, the velocity must vanish 
at the lower limit. 

In this model, and in the one examined next, stability it physically obvious. The virial 
theorem has no interesting application in either case. 



3.2. The polytropic column 

Let us leave the parameter k = p/T'^ free and fix the value of n. This conforms to 
the usual approach when the temperature is fixed by edict, but it is consistent with our 
formulation if n = 3 only, or in the case that the effect of radiation is neglected. We study 
the stability to vertical perturbations. 

The static solution of (2.2) is 



cT = X-gz, c:=n{l + logk). 



A first order perturbation satisfies 

^ + 5X = c5T, thus v=-c5T'. (3.2) 



The equation of continuity gives 

p=-(vpy, p=-{vp)'. 

Let X = X/g — z, < x < Xjg and let /' = dj jdx from now on. Solutions of the type 
(5T = exp(iu;t)/(a;) satisfy the equation 

[x'^bT')' + —{x"bT) = 0, ^ nuJ^jg. (3.3) 

The solution that is regular at the origin of x (the top of the atmosphere) is 

(5T= oFi(n, -iy2a;)e^'^*. 



(Another solution may be obtained by the transformation n — > 2 — n, a; — > 1/x.) The gen- 
eralized hypergeometric function is positive for positive argument and it oscillates around 
zero for negative argument. 
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From (3.3) one obtains 

v'^ j x^'-^STfdz = - J 5T{x'^5T')' = -x^'STST' 



A/a 



+ 



J x^'iST'f. 



Boundary conditions . If we fix 5T = at the bottom of the column the boundary term 
vanishes, provided that the upper hmit makes no contribution, and in this case the mass 
would be preserved. There are solutions for u'^ > only, oscillatory in time. 
The natural boundary condition is the preservation of the mass, thus 

6M = J Spdx = nk j T'^'^STdx = 0. 

This may happen for a discrete set of positive values of u^. For negative values of z/^ the 
integrand is definite so that it can not happen. The calculation is valid only in the case 
77, = 3 (for all n when the effect of radiation is neglected); this atmosphere is stable. 

The problem can be converted to a standard boundary value problem by rescaling of 
the coordinate. 



III. The polytropic gas sphere. The hamiltonian 

Here we study the self gravitating polytropic gas. A correction is needed in the 
expression for the lagrangian, and we need to take care with respect to the definition of 
the gravitational potential. 

The gravitational energy of a system in mutual interaction is the following functional 
of the mass density p, 

EM = -f / d>.d^.:^P^. (3.4) 

It vanishes in the limit that the gas is diffused over an infinite volume and it is negative 
for all other density profiles. 

This potential energy contributes to the hamiltonian so that we have to include a 
term —Eg[p\ in the lagrangian. It contributes to the equation of motion that comes from 
variation of p a term —pcf), as in the other cases examined, but now 

0[p](f) = m = - [ d'x^^ (3.5) 
J \x X \ 

Our text books show how to evaluate the potential in the case of a spherically sym- 
metric distribution, when one writes p{x) p{r), (l>{x) 4>{f)i with 

0(r) = -G j sin 9' dr' d9' dcf)' , z^r^+r''^- 2rr' cos 9', z > 0. 

The trick is to replace (r', 9') by (r',z) as independent variables of intgration, using 

zdz = rr' sin 9' d9'. 
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The integral becomes 

0(r) = -2ttG J p{r')^-dr'dz = -2ttG{I+ + /_), 
where /+ is the contribution from the region r' <r and /_ is the rest. Now 

I+ = - p{r')dr'dz = - / pir'y^dr', 

i-r+r' poo 

I- = - p{r')dr'dz = 2 p{r')r'dr', 

^ J r' — r */ r 

4)'{r) = 477^4 ! p{r')r'^dr'. 



It results from this that 



the mass contained within the distance r from the center is 



M{r) = / p{x')d-'x\ 
I ff.1 



and the result of the calculation is that 

'^W^^' '^(^) = — 

The potential is negative, increasing from the origin to the boundary and on towards zero 
at infinity, since for r > R, (j){r) = —GM./r. 

We shall need 3 important integral forrmulas, the first two are found by by integration 
by parts. 

Eg = -A-K f lp(j)r^dr = \ j M'(j)dr = ^M0(i?) - }- f M(j)'dr, (3.6) 
Jo 2 2 J 2 2 J 

— / M(^ dr = — / — ^dr = - / GMH-) dr = / —MdM. (3.7) 

2j ^ 2 J 2 J 2 R J r ^ ^ 

The third is obtained by Eddington (1926) by using the polytropic relation between 
T and p. We need the equation of motion, 

A-0 = cT, c = 7^(log^ + ^), 

ko 3 

but only to extract the boundary conditions. Since T{R) = it tells us that A = 4>{R)- 
This field, namely 

:= 0(i?) - 0, 
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is Eddington's gravitational potential. The polytropic relation is seen as relating ip and p 
and gives the result, here expressed in terms of 0, is 



/ 



-MdM = — ^ / M(t)'dr'. (3.8) 
r n+1 ' ^ ^ 



Combining (2) and (3) we obtain the data that will allow us to determine the value of the 
hamiltonian, 

jM4>dr=——, J-MdM = — — 

and 

The hamiltonian density is 

h = pvy2 + Ip4> + TZTp log k + 



or, in view of (2.5), 



1 1 3c 

h=-v^ + -p(l>+-Tp, c = 7^(- + log/c), 



whence the hamiltonian (= total energy) is 

H = J d^xh = ^9 + J d^^i^^"^ + ^Tp). (3.10) 
Finally, the equation of motion \ — (j) = cT yields 

cTd^x = AM + 2Ea = — 

5 — n R 

Hence (4.12) reduces, in the static case, to 

3 GM^ 3 n + 1 GM^ 3 n - 3 GM^ 

H = h = . (3.11) 

5-n R 45-n R 45-n R ^ ' 

The conclusion that is drawn from this formula is that the equilibrium configuration is 
stable only if n < 3. 

Remark. It is useful to calculate the relative contribution of radiation to the total energy; 
this can be done exactly when n = 3, using the equation of motion (2.5), 

/, ^ «x 4a 
7eA;log--3 

/Co 3 
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We already gave a numerical example: in the case of the Earthly atmosphere log(A;/A;o) — 
n ^ 10~^°. In general, the internal energy, with radiation included, can be expressed on 
shell as 

U = nRTp + -T^ =: Ugas + Urad- 

and it is natural to ascribe the first term to the gas and the other to radiation. Defining 
the quantity /3 as Eddington does, namely 

/3 — Pgaa / iPgaa Prad) — '^gas / {'^gas ~l" '^rad)-i 

we obtain 

_ log(A:/A:o) - 3 

^ log(/c//co) + 1 ' 

which in the numerical example is approximately .25 x 10""'^°. 

Historical review. 

According to Kelvin (1800), the conclusion that n = 3 is critical for stability was 
announced by Perry, but Perrys paper could not be found. The same result was obtained, 
in the same year, by Ritter (1880) and Betty (1880). Ritter's calculation has often been 
repeated, but the original is by far the most coherent; it is very close to the one that we 
have given above. 

Ritter does not formulate a lagrangian or a hamiltonian dynamics. He defines the 
total energy as the sum of the internal energy and the gravitational energy and postulates 
that this quantity is conserved. He did not have a dynamical theory in which the total 
energy is a well defined functional of the dynamical variables, but it is not a surpris e to 
find that such a dynamical theory exists. 

Later writers are much less careful. Eddington repeats all of Ritter's calculations, 
including the derivation of (3.8) that we have attributed to Eddington because his calcula- 
tion is more easily accessible. But Eddington does not make use of the expression for the 
thermodynamical internal energy; instead he appeals to the interpretation of the gas as 
a collection of particles with mutual gravitational interactions. His claim that the result 
applies when radiation is taken into account is valid only in the case that n = 3, as the 
reader can easily verify. Some modern writers give Ritter's calculation the attention that it 
merits (Weinberg 1972), but some others take shortcuts that make their conclusions far less 
compelling. Kippenenhahn and Weigert (1990) reach the conclusion that n = 3 is critical 
on the basis on an expression that is simply defined to be the total energy. Chandrasekhar 
(1958) also discusses this problem. 

The above calculation made use of the polytropic relation, so it is valid in two cases: 
for all n if the Stefan-Boltzmann term is neglected, and in the case n = 3 whether the 
Stefan-Boltzmann term is included or not. Eddington claims the result for all n, but if 
n 7^ 3, under the condition that the radiation pressure is a fraction of the total pressure, 
fixed throughout the star, which is not verified in our model, and under the additional 
condition that the polytropic relation between temperature and density is valid, which is 
also not verified. 
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All this leaves unsettled the question of stability for the case that the radiation term 
is included. However, if there is a unique, critical value of n then it must be n = 3, since 
Ritter's calculation applies to this special case, even when the energy of radiation is taken 
into account. 



IV. The polytropic gas sphere. Stability 

Here we shall invetigate the stability of the polytropic sphere directly, by perturbation 
theory. 

We use the lagrangian 

£ = p(6 - - 0/2 + A) - UTp log k + ^T^, k := p/T^. (4.1) 

Variation with respect to T gives 

An 

n{\ogk - n) ^ —T^ / p. (4.2) 

With n = 3 this makes k a constant, and log A; = n when radiation is neglected. In the 
remainder of this section, we set, for all values of n, 

p = /cT"^, k constant. 

This is the usual polytropic relation used by Eddington and others, but it is consistent 
with (4.1) only when n = 3. The remaining dynamical equations are 

i)' + cT\ p + r-2(rV)' = 0, 



Dt 

47rGp = r-^{r^(P'y, p = kT''. 



• The static solution. Eliminate (phj cp' = —cT' and change variables, setting r = x/a, 
a constant, Poisson's equation becomes 



where the prime now stands for differentiation with repect to x. Set f{x) = T{x)/T{0) 
and a = \/ AttG / cT(Qi) so that finalUy 

a^'r + (^Y)' = 0, /(o) = i, /'(o) = o. 

The solution decreases monotoneously to zero at x = X, this point taken to be the surface 
of the star. At the outer limit f{x) (x X/ x — 1 + o{X — x)"^ . The integration is done easily 
and accurately by Mathematica, especially so for integer values of n. The radii are, for 
n = 2:X = 4.355, n = 3 : X = 6.89685636197, n = A:X = 14.9715. 
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• For the fluctuations we assume harmonic time dependence, then the equations are 

-uj\'^5p = {r^p {S(f)' + cST'))', Sp = nkT'^-^ST, (4.3) 

AirGr^Sp = {r'^64>'y. (4.4) 
Introduce the function 6Ai — r'^64>' . Eq.s(4.3-4) then take the form 

-——SM = pdM + r pc5T + constant, 

where the constant can only be zero, and 

(47rG)r2(nA;T"-^5T) = 5M' , 

Ehmination of 5T leads to 

Changing the scale as before we get 

-.^6M = nM + l^y-i^)', = sg^. (4.5) 

The crucial point is the choice of the correct boundary conditions, at x = as well as 
the outer surface {x — X). At the center the solutions take one of two forms, 1 + Cx'^ + 
which is unphysical, or else + Cx^ + ... Accordingly we set 

SM{x) = x^g{x), g{0) = l, ^'(0) = 0. (4.6) 

The boundary conditions at the outer boundary are determined by the fact that the mass 
is conserved, 

SM = SM{X) = 0. 

The equations then imply that the zero is of order n. With these boundary conditions (4.5) 
becomes a well defined Sturm-Liouville problem with an essentially self adjoint, second 
order differential operator. 

Numerical calculations with the help of Mathematica are not difficult in the case of 
integer values of n. It is found that, when n — 2 and for n — 3, dA4{X) is positive in 
the whole range, for all negative values of and for positive values below a limit Vq that 
is about .06 for n = 2 and compatible with for n = 3. The latter is the first, nodeless 
solution of a sequence of solutions that we have not determined in detail. The function 
falls to zero at the surface, where there is an nth order zero. Above this lowest value of 
is a discrete set of other values of u"^ at which the boundary condition is satisfied. 

At the special value n = 3 the 'ground state', the lowest value of u"^, has approached 
very close to zero. 
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Polytropes with n = 4 are unstable. We have searched for harmonic solutions with 
negative values of . The value n = 4 was chosen because it is the only integer in 
the interesting range, and because Mathematica is much more managable in this case. 
(Accuracy is lost when non integral powers of negative numbers appear at the end point.) 
There seems to be a discrete, decaying nodeless mode with = —.015796, but a bifurcation 
at this point in parameter space makes the conclusion uncertain. We carried the calculation 
to 15 significant figures in but solutions do not converge towards a function that vanishes 
at the surface. To overcome this difficulty we reformulated the problem in terms of the 
variational calculus. The "solution" found for — —.015796, truncated near both ends, 
was used as a trial function, to show conclusively that the spectrum of extends this far. 
Among many papers on this topic we mention Cowling (1936) and Ledoux (1941). 



V. The case n = 3 

This case marks the boundary between stable and unstable polytropes. The equations 
are conformally invariant and a time independent solution is found by an infinitesimal 
conformal (homology) transformation, 

Sf = rf' + f. (5.1) 

This does not represent an instability, but a "flat direction" , a perturbation from which the 
system does not spring back, nor does it run away. There must also be a second solution, 
linear in t, of the form 

Sf = t{rfy, Sp = t{rp' + 3p). 

The equation of continuity becomes rp' + 3p + vp' + r~'^{r'^v)' p = 0, whence v = —r. 

This linear perturbation is the first order approximation to the exact solution found 
by Goldreich and Weber (1980), of the form 

/(^' = ■^(^)' ^ = r/a{t). 
a{t) 

The continuity equation is solved hy v = ax; thus $ = — (d/a)(r^/2), and 

- = -aax'^/2 = cT + (j). 

This leads to ^ 

4> a{t)4> oc / + «;a^ax^/6, k = 3k^^^/c, 

and Poisson's equation becomes 

~ 1 ~ —K 

f + -^ix^f')' = ^a^axV6 = -na^'d = A, constant. (5.2) 
a;^ x^ 

There is a first integral, 

—b? — \/a = constant. 



17 



Rescaling of t and a reduces this to one of three cases 



d = a/1 + 1/a, a = a/1 — 1/a, d = 1 / ^a, 

but only the first is compatible with analyticity at t = 0, thus 

t = ^/a^/l + a — arcsinh^a. 

Setting a = 1 + 6 we find 

t = ^/T/2{b - /2) + oib^) 

The factor a{t) is zero at a finite, negative value of t and increases monotoneously to 
infinity, passing through 1 at t = 0. We can of course reverse the direction of fiow of t to 
get collapse in the finite future. 

Eq.(4.5) was solved numerically (Goldreich and Weber, 1980). The solution is similar 
to the solution of Emden's equation, just prolonged a little at the outer end, so long as 
< A < .00654376. For larger values of A the distribution does not reach zero and increases 
for large r. For simillar studies of collapsing, isothermal spheres see Hunter (1977) and 
references therein. 

It is confirmed, therefore, that the polytrope with tt, = 3 is not stable. Suitably 
perturbed, the star may expand or collapse, until the higher or lower density causes a 
change in the equation of state. Among may papers on collapse we may mention Arnett 
(1977) and Van Riper (1978). 
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